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Abstract. This article is an attempt to provide a link between the quantum 
nonequilibrium dynamics of cold gases and fifty years of progress in the low- 
^ ' dimensional quantum chaos. We identify two atomic systems lying on the 

C^ ' interface: two interacting atoms in a harmonic multimode waveguide and an 

3 ' interacting two-component Bose-Bose mixture in a double-well potential. In 

^y^, particular, we study the level spacing distribution, the wavefunction statistics, the 

eigenstate thermalization, and the ability to thermalize in a relaxation process as 
C^ ■ such. 



PACS numbers: 67.85.-d, 37.10.Gh, 05.45.Mt, 05.70.Ln 



Two simple systems with cold atoms: quantum, chaos tests and nonequilibrium dynamics2 

1. Introduction 

In this article, we identify two simple atomic system lying on the interface between 
the nonequilibrium dynamics of cold gases and low-dimensional quantum chaos. The 
system of two short-range-interacting atoms in a harmonic waveguide is compared to 
the well studied Seba billiard Ijy. A two-component interacting Bose-Bose mixture 
in a double-well potential provides a realization of a quantum particle moving on a 
four-dimensional surface of a tensor product of two spheres under a non-integrable 
perturbation. 

The waveguide problem allows for an exact analytic solution, in spite of the 
absence of a complete set of integrals of motion. The double-well problem reduces to 
a numerical diagonalization of a small, a-few-hundred-by-a-few-hundred matrix, a few 
minute calculation on a moderate laptop. 

To assess the degree of quantum chaos, we utilize two standard quantum- 
chaotic measures: level spacing statistics [21 [3] and the statistics of the values of 
the wavefunction in the coordinate representation [H [S] . The former predicts that the 
spectra of quantum-chaotic systems have the statistical properties of the spectra of 
random matrices. The latter dictates the Gaussian statistics for the wavefunction. 

From the point of view of the nonequilibrium dynamics, we study both the degree 
of the eigenstate thermalization and the actual thermalization in an expansion from 
a class of realistic initial states. The eigenstate thermalization [51 [71 |H1 [SI [ID] - the 
suppression of the eigenstate-by-eigenstate variance of quantum expectation values of 
simple observables - provides an ultimate upper bound for a possible deviation of the 
relaxed value of an observable from its thermodynamical expectation, for any initial 
state in principle. However, recently a new direction of research has emerged: quantum 
quench in many-body interacting systems pT | fT2 l [T3 l fT4 l [T5 l [T6 l [T7 l [T8 | fT9 l [20 l [2T ] . 
In this class of problems, the initial state is inevitably decomposed into a large 
superposition of the eigenstates of the Hamiltonian governing the dynamics, and the 
discrepancy between the result of the relaxation and thermal values is expected to be 
diminished. 

Note that a related comprehensive study of the relationship between the quantum 
chaos and thermalization for atoms in optical lattices has been performed in Ref . |22j . 

2. Tv^fo bosons in a circular, transversely harmonic multimode waveguide 

In the case of the relative motion of two short-range-interacting atoms in a circular, 
transversely harmonic waveguide [53] , the full Hamiltonian 

H = HwG + Vfh (1) 

can be split onto an integrable, "free motion" part and a non-integrable perturbation. 
The unperturbed Hamiltonian is given by the sum of the longitudinal and transverse 
kinetic energies and the transverse trapping energy U = /ia;^p^/2, 

where uj± is the transverse frequency, Ap is the transverse two-dimensional Laplacian, 
12 = to/2 is the reduced mass, and to, is the atomic mass. We will assume periodic 
boundary conditions along z, with a period L. In what follows, we will restrict the 
Hilbert space to the states of zero z-component of the angular momentum and even 



HwG^- — 7n^-—^P + U-fi^± , (2) 
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under the z ^^ —z reflection. Note that the zero-range interaction has no effect on 
the rest of the Hilbert space. The non-interacting eigenstates are products of the 
transverse two-dimensional zero-angular-momentum harmonic wavefunctions, labeled 
by the quantum number n > 0, and the symmetric plane waves cos(27rZz/L), / > 0. 
The unperturbed spectrum is therefore given by Eni = 2hiu±n + fi,^(27r^/L)^/(2/i) 

The interaction Vfh in ([2|) couples the transverse and longitudinal degrees of 
freedom. We assume the interaction to be of the Fermi-Huang type. It can be 
symbolically represented as a separable rank I perturbation |24] : 

VFH^v\c){ni (3) 

with the formfactors |£) = S3{r) and {TZ\ = 63{r){d/dr){r ■). Here, the interaction 
strength is 1^ = (2iTh^as/ p,), where Og is the three-dimensional s-wave scattering 
length. The rigorous definition of the operator (j31) reads 

(x|VFi/|V') = ^X*(r = 0)lim^(r^(r, 0, ^)) . (4) 

r— i>0 or 

Following the the Refs. [Ml [T], it can be shown that the operator ^ produces a 
member of a (parametrized by V) family of self-adjoint extensions of the "free motion" 
Hamiltonian ([2|). 

Note that the model H]), with an unbounded along z motion, was first used to 
derive the effective one-dimensional atom-atom interaction in a monomode atom guide 
pS] . It was further extended to the multimode regime in Ref. [55] , where the intermode 
kinetic coefficients were derived. In this article, we consider a bound-motion analogue 
of the multimode guide of Ref. [2^, suitable for a study of the relaxation dynamics 
and chaotic properties. 

Even though the perturbation of the form ([3]) destroys the complete set of integrals 
of motion of the unperturbed system, the problem can be solved exactly. The solution 
for the full system can be obtained as follows. The action of the perturbation on each 
wavefunction produces the left formfactor \C). As a result, any eigenstate |a) of the 
perturbed system functionally coincides with the Green function of the unperturbed 
Hamiltonian taken at the energy of the eigenstate Ea. The following eigenenergy 
equation then applies, along with an expression for the eigenstates: 

n n 

where \n) stands for an eigenstate of the unperturbed system of energy £',-j. Similar 
expressions were obtained in Refs. [TJ [37] for the case of the Seba billiard and 
its generalizations. Solutions of this type, for the integrable case of a spherically 
symmetric harmonic confinement, are presented in Ref. [28] . 

Substituting the above non-interacting eigenstates and eigenenergies to the Eq. 
([5]) and summing over /, one obtains the following expression for the interacting 
eigenstates: 

,7^^ V £« - An sm Vca - An 

where the rescaled energy Cq. is given by £„ = \E/{2fiw±), A = (L/a±) is the 
aspect ratio, C = z/L— 1/2, S, = ip/a±) , Ln {£,) are the Legendre polynomials, 

1 /9 

and a± = {h/fiLj±) is the size of the transverse ground state. The normahzation 
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constant C^ is expressed as 

/ rinf A / c — An 1 \ 

(7) 



Cq 



a^y/irL 



E 

n=0 



(eq — Aji) 3/2 (f^ — An) sin^ V^a — ^' 



Extracting the regular part of the wavefunction ^ via the procedure given in 
Refs. [3Sl[3n]j we arrive at the foUowing transcendental equation for the eigenenergies 



/-r v^ cot Veg - An + t /I ea\ _ a± 

;^ y/ea - An V2' A/ a^ 

where C(j^, a;) is the Hurwitz (^-function (see [29]). The summands in the sums Eqs. 
©j ^, and (j8]) decay exponentially with n, leading to the fast converging series. 
Note that the imaginary parts of the two terms in the left hand side of the ([S]) cancel 
each other automatically. Similar solutions were obtained for two atoms with a zero- 
range interaction in a cylindrically-symmetric harmonic potential [30' (that system 
was analyzed numerically in Ref. |31j). Higher partial wave scatterers were analyzed 
in the Ref. [32] . 

At rational values of A/tt^ the unperturbed energy spectrum shows degeneracies 
that are not fully lifted in the full [deduced from (|S])] spectrum. To minimize the 
effect of the degeneracies, we, following Ref. [1], fix the length of the guide by 
{L/a±)'^ = A = 2(j)grTr'^ ~ 9774, where cpgr = (1 + \/5)/2 is the golden ratio. 

3. A tv^ro-component interacting Bose gas in tv^ro coupled potential wells 

Next, we are going to consider a two-component Bose-Bose mixture trapped in a 
two- well potential. 

The second-quantized Hamiltonian is 

H = Ho + V (9) 

^0 = 2 E E (t^yi^iSxC^LSL + U^jBclB.CJRC^RCjRj (10) 

i—a, b j—a, b 

f = -2 J ^ (cLc,« + cl^qj , (11) 

i—a, b 

where cjj^(c|^) is the creation operator for the atom of the i-th type in the left(right) 
well, a and b denote the two species, and the factor of —2 in the definition of V is used 
to ensure the consistency with the conventional definition of the hopping constant in 
the theory of multi-well lattices with periodic boundary conditions. 

We fix the coupling constants f/ijL(i?) to UaaB. = 0.3041 [/, UabR = 
lM10U,UbbR = 0.3684 [/, C/aaL = 0.3085 [/, C/a6L = lA574:U,UbbL = 0.4524 C/, 
where U is the overall coupling strength. The relative strengths of the interaction 
were chosen once and for all from a uniform random distribution between and 2. 
The same set was used in all further numerical experiments. 

The set of deviations of the number of atoms of a given type in the right well 
from the half of the total number of atoms of this type constitutes a convenient set 
of quantum numbers: Sa = NaR - Na/2, s^ = NbR - Nb/2, where NaR = c\^c^^ 
and NbR = c\jj^c^^. The observable Sa has been chosen to be the principal observable 
of interest in the eigenstate thermalization and relaxation studies. Below, we fix the 
numbers of atoms to Na = Nb = N = 20. 
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Analogously to the one-specie case (33] , the Hamiltonian ^ is integrable in both 
J :$> U and J <t:U limits. 

In the J <^ U limit, the eigenstates are represented by the Fock states 
(eigenstates of Hq, see PTI)1 ) where number of bosons of each type in each well is 
fixed. For our choice of the coupling constants, the energy as a function of (sq, Sf,) 
constitutes a hyperboloid. At low energies, the spectrum of Hq is dominated by 
the eigenstates where the atoms of different type are localized in the opposite wells: 
(sq, Sb) ~ (±A^, TN). At high energies, the atoms prefer to localize in the same well: 
{Sa, Sb) « i±N, ±N). 

In the opposite limit J ^ U, the good quantum numbers are given by the numbers 
of particles of each type in the even and odd (with respect to the wells) one-body states. 
In this limit, the mean- field Hamiltonian 

i?,nf =V+ (Ho) (12) 

N 

= V+-{{N- 2){UaaL + UaaR + UbbL + Ubbn) 
o 

+ 2N{UabL + UabR)) , 

with the hopping as the principal component, generates a good approximation to the 
energy spectrum. The spectra of the full Hamiltonian H ([9]), the pure-interaction 
Hamiltonian Hq (jlOl) . and the pure- hopping- mean- field Hamiltonian H^f are shown 
at Fig. inh- Note the degeneracies in the spectrum of the pure-hopping Hamiltonian 
Hmf- Recall that in our model, the hopping constants are the same for both types of 
atoms; as a result, each energy level corresponding to the total of N+ atoms in the 
even state {1/V2){\L) -\- \R)) is {N - \N+ - N\ + l)-fold degenerate. 

The point a ~ 190 where the exact, pure-interaction, and pure- hopping 
eigenenergy curves coincide corresponds to the maximal depolarization, {sa(b)) ~ 
^ (b)L^a(b)B) ^ ^' ^^^^re neither of the two integrable limiting models apply. This is 
the point where the motion is expected to be "maximally chaotic" . 

The attractive feature of the two-well system is that it can be mapped to a 
system consisting of a single particle moving on the four-dimensional surface of a 
tensor product of the two two-dimensional spheres. The map is |sa, Sb) — > \{La — 
Na/2, ma = Sa), {Lb = Nb/2, mb = Sb)), where La{Lb) is the total angular momentum 
on the first(second) sphere, and ma{mb) is the corresponding projection to the z-axis. 
The corresponding Hamiltonian reads: 

H = ]^{UaaLLA^-.a ^ 1) + 2UabLL,aL,b + ■ ■ ■) (13) 

i—a, b 

The total angular momenta (Ca)'^ = La{La + 1) and {Cb)'^ = Lb{Lb + 1) commute 
with the Hamiltonian (J13p . In the many-boson picture, these momenta are related to 
the number of bosons of each type as ia(6) = A^o(&)/2 (both assumed to be even). For 
the one-specie case, similar mapping was employed in Ref. |33| . 

From the single-particle point of view, the Hamiltonian (J13p describes an external 
perturbation of the free motion along the bi-sphere. For clarity, we do not include the 

constant "free" kinetic energy ex {CaY + i^b)'^ in the Hamiltonian (fT3|) . 
In the coordinate representation, a given state 

\^) = y] ('5'^' ■Sf'IV')|Sa, Sb) 
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Figure 1. Two bosons in a multimode guide. Nearest-neighbor spacing 
distribution averaged over 7 X 10® eigenstates in the globally unitary regime 
as/aj_ = 10® and in the regime of no interactions. The Gaussian Orthogonal 
Ensemble prediction, the Seba distribution, and the Poisson distribution are 
shown for comparison. Here and below, {L/aj_) = 2<f>gr'n'^ ~ 9774. 



maps to a two-particle wavefunction (defined at the points of a tensor product of two 
unit spheres) given by 

V'(0a, 4>a, Ob, 4>b) = (14) 

where Yi^mid, <p) are spherical harmonics. 



4. Standard quantum-chaotic tests 

^.1. Level spacing statistics 

The commonly-used measure to access the degree of quantum chaos the shape of 
level spacing distribution [3]. As usual, the starting point of the procedure is to 
determine the mean level density to be subsequently removed from the spectral data. 
The resulting energy spectrum — the so-called unfolded spectrum — does not contain 
any information specific to the system besides the dimensionless, homogeneous across 
the spectrum level spacing statistics. This statistics is the principal object of the 
investigation. 

For the case of the waveguide, the mean level density was obtained from the 
semiclassical approximation to the unperturbed spectrum of the system. Such a 
simplification was possible due to the fact that for a singular perturbation, the density 
of states is exactly the same for the interacting and non-interacting systems. For the 
two-well system, the mean level density was approximated by a fit of the true energy 
spectrum to a seventh degree polynomial. 

The integrable systems exhibit the Poisson statistics of the level spacings, with 
a peak at zero and an exponential decay at large spacings. One would expect to this 
type of statistics in a model the energy levels appear under no constraint, besides the 
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N3=Nb=N=20 I a= [150:250] 



J/U = 2. 

GOE 

J/U = 

Poisson 







0.5 1 1.5 2 2.5 3 3.5 
Unfolded level spacing, Ae 



Figure 2. Forty bosons in two wells. Nearest-neighbor spacing distribution 
averaged over 100 eigenstates in the "maximally chaotic" region 150 < a < 250. 
Number of atoms of both type is Na = N^, = N = 20. The J/U = 2 point 
corresponds to the empirically broadest fraction of the "chaotic" eigenstates. The 
case of no hopping is shown as well, along with the Gaussian Orthogonal Ensemble 
and Poisson predictions. Here and below, the individual coupling constants are 

UaaR = 0.3041 t/, UabR = 1.04101/, UbbR = 0.3684 t/, UaaL = 0.3085 t/, UabL = 

1.4574(7, UbbL =0.4524(7. 



given a priori mean level spacing. Indeed, in integrable systems, the energy levels 
corresponding to different sets of (additional to the energy) integrals of motion can 
approach each other infinitely closely. It turns out that this property along is sufficient 
to ensure the Poisson statistics [21 [3] . 

On the contrary, in a quantum system with no non-trivial integrals of motion, the 
energy levels tend to repel. This repulsion leads to the appearance of a hole, around 
zero spacing point, in the level spacing distribution. It is remarkable however that 
generally, quantum chaotic systems with real Hamiltonians exhibit exactly the same 
level spacing distribution, namely the one of the Wigner-Dyson type. An ensemble 
of orthogonal matrices with a Gaussian distribution of the matrix elements (the so- 
called Gaussian Orthogonal Ensemble (GOE)) serves as a paradigm for this class of 
Hamiltonians [21 [3] ■ 

A multi-dimensional harmonic oscillator constitutes a remarkable exception from 
the integrable vs. chaotic dichotomy [34J. The unfolded distribution function in this 
case vanishes at the zero spacing and shows a sharp decay at large spacings. Indeed, 
we have found that our waveguide system, with a harmonic confinement transversally 
and a periodic-boundary-conditions box longitudinally, is more sensitive to the effect 
of a non-integrable perturbation than a seemingly more elegant elongated three- 
dimensional harmonic oscillator; thus our preference for the former. 

The waveguide system shows systematic deviations from the GOE Random 
Matrix Theory predictions [2 . Note that the unattainability of a complete quantum 
chaos in singular billiards has been already addressed in the context of the Seba billiard 
- flat two-dimensional rectangular billiard with a zero-range scatterer in the middle [1] . 
Even though the Seba billiard shows some signatures of the quantum-chaotic behavior, 
both the level statistics |351 127| I36j and the momentum distributions in individual 
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Figure 3. Two bosons in a multimode guide. The probability distribution for 
value of the eigenstate wavefunction in the region close to the waveguide axis. 
The eigenstate number 67600 corresponding to the energy Ea ~ 4.24 X lO^fitJx 
is shown. Since the motion is restricted to the cylindrically symmetric states 
the reduced wavefunction V(Pi ^) = (l/\/2''rp)V'(Pi ^) wa-s used. The Gaussian 
distribution corresponds to the Berry conjecture. The interaction strengths 
as/aj_ = 10^, lO"'^, and 10~^ correspond to the unitary, intermediate, and non- 
interacting regimes respectively. The rest of the parameters is the same as at Fig. 

m 



eigenstates |37] exhibit substantial deviations from the conventional quantum- chaotic 
behavior. Our waveguide system belongs to the same class as the Seba billiard: it is 
represented by an integrable system perturbed by a separable rank-one perturbation. 
For both systems, the eigenstates have the "Green's function" form ([Sj- 

The results of the study of the level spacing distribution in the waveguide system 
are fully consistent with the analogous results for the Seba billiard [351 HZ]- One 
can show that for Os 3> aj^ (globally unitary regime), the scattering strength reaches 
the unitary limit for all eigenenergies of the full Hamiltonian. In this regime the 
distribution quickly converges to the Seba distribution [35] • This distribution does 
show a gap at small level spacings but fails to reproduce the Gaussian tail predicted 
by the GOE. At small a^ the level spacing distribution tends to the Poisson one (see 



Fig. [Ij. For the energy range E > 10Qhuj± and the aspect ratio {L/a±) 



24>g 



the Seba distribution is approached at Us > 10~ a±. 

The level spacing statistics for the two-well system is consistent with the GOE 
prediction (Fig. [2]), even though the statistics is very poor. 



4-2. Wavefunction statistics 

Another important measure of quantum chaos is the statistics of the wavefunction. 
Divide the coordinate space of a quantum-chaotic system into regions larger that the 
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Figure 4. Forty bosons in two wells, (a) The probability distribution for 
real part of the eigenstate wavefunction on the whole surface of the bi-sphere. 
The distributions for the individual eigenstates were further averaged over an 
150 < a < 250 window. The distribution is heavily peaked in both integrable 
limits, J <^ U and J ^ U. (b) Entropy of the distribution as a function of the 
hopping constant. The entropy converges to the Berry conjecture prediction in 
the region between J/U ~ 1 and J/U ~ 5. For both (a) and (b) the rest of the 
parameters is the same as for the Fig. [2] 
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Figure 5. Forty bosons in two wells. The real and imaginary parts of the four- 
dimensional eigenstate wavefunction living on the surface of a tensor product of 
two unit spheres. The "maximally chaotic" state a = 190 is used. The hopping is 
fixed to J/U = 2. Even though in both J <^U and J ^ U limits the eigenstates 
are strongly polarized, for J/U = 2 no preferential direction is visible. The rest 
of the parameters is the same as for the Fig.|2] To enhance contrast, the shading 
density is made to be proportional to the forth power of the real and imaginary 
parts respectively. 



de-Broglie wavelength but small enough so that the spatial density does not change 
substantially across the region. According to the Berry conjecture 01 [5], the coordinate 
wavefunction should behave, within each region, as a sample from an ensemble of large 
superpositions of plane waves with random coefficients. 

For the waveguide system, we verified that in the globally unitary regime. 
Us ^ a±, the statistics of the point- by-point variations of the eigenstate wavefunction 
in the spatial representation is indeed close to the Gaussian one (Fig. [3]), according to 
Berry's prediction and a particular observation in the case of the Seba billiard 1 . 

The two-well system also behaves according to the Berry conjecture. Fig. 0^. 
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Figure 6. The same as on Fig.[5]but in the regular regime, J/U = .2 . 



presents the distribution of the real part of the wavefunction (J14p taken at random 
points on the bi-sphere and averaged over the eigenstates between a = 150 and 
a = 250. For J = 2 the distribution is indistinct from the Gaussian. Fig. lib represents 
the dependence of the entropy of the distribution on the hopping constant. The 
entropy converges to the Berry conjecture prediction in the region between J/U ~ 1 
and J/U « 5. 

Fig. [5] shows the wavefunction of the eigenstate a — 190 at hopping of J = 2. 
This state corresponds to the crossing of all three spectra of the Fig. [9^ where the 

expectation values of the projection of the angular momentum Ca (along with Cb) 
on both z and x axes is small. Accordingly, the nodal structure of the wavefunction 
appears completely random, with no preferred direction. For comparison, Fig.|6]shows 
the same state but for weak hopping (J = .2) where the motion is expected to be 
regular. 
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Figure 7. Two bosons in a multimodc guide, (a) Quantum expectation 
value of the trapping energy U for every hundreds eigcnstate of the 
unperturbed Hamiltonian |(2]| (blue circles) and the full Hamiltonian ((ij (grey 
squares), (b) The time average of the transverse trapping energy C].gja,x = 
li'^tmaj;— >oo(l/*™.a2:) /o""'^{'!/'(*)|f?IV'(*)) ^ftcr rclaxatiou from an initial state as 
a function of the energy of the initial state E [see Hamiltonian jT}] . Four families 
of the initial states are considered. Long-dashed (green) line: state JTSj with 
no = 0, 5 = 0.99, and the energy being controlled via scanning Zq; short-dashed 
(blue) line: the same as the previous one, except for 5 = .1; dotted (purple) line: 
the state I I15I I with <5 = 0.1 and Zq = and the energy control via no; solid (red) 
line: state I I16II with K2 = 2ki, 5 = 0.99, Iq = and the energy control via fti. 
All the values of the trapping energy are shown with respect to its microcanonical 
expectation value Utf^f,,.^ ~ E/3. The rest of the parameters is the same as at 
Fig.IIl 
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Figure 8. Two bosons in a multimodc guide. Probability distribution of the 
quantum expectation values of the trapping energy U for individual eigenstates 
of the full Hamiltonian JlJ in the globally unitary regime (solid line) and those 
of the unperturbed Hamiltonian l(2} (dashed line). The rest of the parameters is 
the same as at Fig. [l] 
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Figure 9. Forty bosons in two wells, (a) Energy spectrum of the Hamiltonian 
(|9]l in the absence of hopping {J/U = 0) and for the hopping J/U = 2. The 
spectrum of the mean-field Hamiltonian H^^f = V + {Hq) (see II12I I. with the 
hopping as the primary Hamiltonian, is shown as well. The point a S! 190 where 
all three eigenenergy curves coincide corresponds to the maximal depolarization 
where the motion is "maximally chaotic". Here, the density is distributed 
uniformly over the bi-sphere with no preferred orientation. The labels a and qq 
correspond to the eigenstates of the full Hamiltonian ^ and the pure-interaction 
Hamiltonian (J/U = limit) respectively, (b) Fluctuations of the observable 
Sa = c],dC^„ — A^a/S. Green crosses correspond to the microcanonical fluctuations 
in the alasence of hopping. Blue crosses correspond to the eigenstate-by-eigenstate 
fluctuations of the quantum expectation values, suppressed according to the 
eigenstate thermalization hypothesis. Red crosses give the infinite time average 
of the quantum expectation for a propagation from an initial state corresponding 
to one of the eigenstates of the Hamiltonian Hq. The fiuctuations of the latter 
are further suppressed due to an additional averaging of the initial state over the 
eigenstates. The microcanonical average is shown as the blue dashed line. 
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5. Eigenstate thermalization phenomenon and the ability to thermalize in 
a relzixation from an initial state 

According to the Eigenstate Thermalization Hypothesis O [7j [U [10], the ability of 
an isolated quantum system to thermalize follows from suppression of the eigenstate- 
by-eigenstate fluctuations of the quantum expectation values of relevant observables 
(a|A|a) over the eigenstates \a) of the system. Indeed, if one assumes that an 
observable A reaches its thermal value in a time evolution from any initial state with 
an energy close to a given energy _E, including the exact eigenstates of the Hamiltonian, 
then it immediately follows that the matrix elements (ai|A|ai) and (a2|^|ft2) must 
be close to each other if the corresponding energies Ei and £'2 are close. 

Consider first the system of two particles in a waveguide. As the principal 
observable of interest, we choose the transverse trapping energy U (see ([2])). In our 
analysis, the potential energy U plays a role of a generic observable in an integrable 
system: it acts on only one of the degrees of freedom and it has both diagonal and off- 
diagonal non-zero matrix elements between the eigenstates of the integrable system. 

Scattered grey squares at the Fig. [7^ show the quantum expectation values of the 
transverse trapping energy, (a|f7|a), in every hundreds eigenstate ([6]) of our system 
(dl, relative to its microcanonical expectation value U therm. ~ E/3 (in accordance 
with the energy equipartition between the two harmonic and one free-space degrees 
of freedom). 

Even for the energies much larger than any conceivable energy scale of the system, 
the variation of (a|[/|Q;) is only marginally lower than the microcanonical fluctuations 
in the unperturbed system ([2]), represented by the blue circles. Figure [8] demonstrates 
that even though the destruction of the integrals of motion in the interacting system 
does narrow the distribution of the diagonal matrix elements (a|C/|Q;) and it does 
shift the peak towards the microcanonical prediction Utherm., the distribution remains 
relatively broad. We trace this effect to the predicted unattainability of the complete 
quantum chaos in the singular billiards |35j |36l |38] . 

Here and below, we assume the globally unitary value of the scattering length, 
a^ = lO^a^. 

On the contrary, the behavior of the two-well system is fully consistent with the 
Eigenstate Thermalization Hypothesis j6l [71 HI [10] (see Fig.[9jD). Here, the observable 
Sa (the deviation number of the type-a particles in the right well from N/2) is the 
primary observable of interest. In the window 150 ^ a < 250, all quantum expectation 
values of Sa in the eigenstates of the full Hamiltonian (|9]) are substantially closer to 
the thermal prediction than their unperturbed counterparts given by the eigenstates 
l^o) — |saOj Sfco) of the Hamiltonian ([TU|. 

Now we are going to address directly the ability of our system to thermalize from 
an initial state \ijj(t = 0)). In this case, the infinite time average of the quantum 
expectation value of an observable A will be given by 

Aclax. = lim {l/t,nax) / {tp{t)\A\tP{t)) 







= ^|(a|V^(t = 0))p(a|i|a) , 

a 

where \a) are the eigenstates of the system (see [TO]). 

First, let us attract our attention to the waveguide system. 
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Consider the following initial state: 

{zm^O)) = Caxcos^e (^ - \C\\ cos(27r/oC) |"o) (15) 



C^ = L-'/' 



1/2 

S 1 vr^ sin (27r/o^) 

4 8^ n^ - {2TTloSf 



Longitudinally, the state is represented by the ground state of a length L5 hard-wall 
box split initially by an ideal beamsplitter with momenta ±2ttIo/ L; the box is centered 
at the maximal interatomic distance. The state is distributed among approximately 
tt/S axial modes localized about I = ±Zo- The initial transverse state is limited to a 
single mode uq. Note that for hq = and 5 « 1 the state is conceptually similar to 
the initial state used in the equilibration experiments |19) . 

The transverse trapping energy U (see ^) remains very far from the equilibrium 
predictions after the relaxation. Fig. [7|d shows the equilibrium value of the transverse 
trapping energy U for three sequences of the initial states of the type P^ . For the 
first two, the energy is controlled via the kinetic energy given by the beamsplitter, 
while the transverse state is fixed to the ground state. Both show an approximately 
-40% deviation of the relaxed value of U from the thermal prediction. Note that the 
initial states of the second sequence have a much greater spread over the longitudinal 
modes than the first one; consistently, the energy-to-energy variation of the relaxed 
values of U is less than for the first sequence. In the third sequence, there is no 
beamsplitting, and the energy is controlled via the initial transverse energy. In this 
case, the deviation from the thermal prediction ranges between -1-50% and -1-5% but 
never reaches zero. 

Another type of the initial state 

(p, zm = Q)) = Caxcos^e (^ - \C\] cos{2ttIoC) 

7r(Kl-K2) J 

allows for an additional spread over the transverse modes. The transverse 
wavefunction vanishes at the waveguide axis. The mode occupation has a minimum 
at n = 0, and the occupation of the ground transverse mode tends to zero when 
Ki < K2 -^ 1. The results for a sequence similar to the third sequence described 
above are shown at Fig. [Tb as well. In spite of the significant differences in the initial 
distribution over both the transverse and the longitudinal modes, the deviations from 
the equilibrium are close to the ones for the third sequence, even though the energy- 
to-energy variations for the former are less than for the latter. 

For our second system, we use the eigenstates of the unperturbed Hamiltonian 
(jlOp as the initial states of the relaxation process. The results are shown at Fig. [3]d. 
Here, each red cross shows the infinite time average of the quantum expectation value 
of the observable Sa (the deviation of the type-a occupation of the right well from 
Na/2) for a relaxation process originating from the state |ao) = \saO, sto}- In the 
window 150 ^ a < 250 the observable thermalizes fully. Recall, that the eigenstate 
themalization occurs in the same window. 

It is instructive to analyze our results from the point of view of the macroscopic 
quantum self-trapping, predicted in Refs. [33l [39l [40] and experimentally observed in 
Ref. [H] , for the single-specie two- well bosonic system. From the classical field theory 
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point of view it is a one-dimensional Hamiltonian system similar to a pendulum. Here, 
the population imbalance between the two wells plays a role of momentum, while the 
relative phase between the wells plays a role of a coordinate. Similarly to the case of 
the conventional pendulum, the initial states of sufficiently high momentum preserve 
the sign of the latter; the time evolution that starts from the states with even higher 
momentum shows the momentum that almost does not change in time. The states 
with the opposite sign of the momentum become mutually inaccessible. In terms of 
the original two-well system, the initial states with the population imbalance higher 
than a certain critical value preserve the imbalance over time. 

The energy landscape for out two-specie two-well system is more complicated. 
Both the cosine-like dependence of the hopping energy on the relative phases and the 
hyperbolic shape of the interaction energy surface as a function of the two population 
imbalances contribute to the potential multi-connectedness of the phase space. At the 
same time, the initial states used to produce the Fig. ^p are the exact eigenstates of 
both population imbalances; they leave the relative phases completely uncertain. 

We have identified a set of absolutely localized pairs of population imbalances, i.e. 
the pairs for which no set of imbalances of opposite signature is classically accessible 
for any pair of initial phases. In particular, the analysis shows that all the unperturbed 
eigenstates with indices ag < 39 or ao > 339 constitute the absolutely localized initial 
states. This result is consistent with the fully quantum calculation that shows a strong 
memory of the initial conditions (see Fig. [9}d). 

6. Summary of results 

In this article, we analyze the behavior of two simple atomic systems within two 
complementary frameworks: quantum chaos and nonequilibrium dynamics. The first 
system is represented by two short-range-interacting bosonic atoms in a circular, 
transversely harmonic multimode waveguide. The second is the interacting two- 
component Bose-Bose mixture in two coupled potential wells. For both systems, 
we study the level spacing statistics, the statistics of the values of the wavefunction 
in coordinate representation, degree of the eigenstate thermalization, and the 
thermalizability in a relaxation from an exited initial state. 

For the waveguide system, the wavefunction statistics is fully consistent with 
the Berry conjecture [1] that predicts a Gaussian distribution for the quantum-chaotic 
systems. However, the rest of the tests show an incomplete chaotization. We trace this 
effect to the previously predicted incomplete quantum chaos in the singular billiards 

[sniiaiMiisz]. 

For the double- well system, we identify a quantum-chaotic region of the spectrum 
that spans, for our set of parameters, about a quarter of the full spectrum. In this 
part of the spectrum, the system exhibits both strong eigenstate thermalization and 
the ability to fully thcrmalizc from an initial state. 
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